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(Abstract) This paper presents a method of solving three-phase probabilistic load flow (PLF) for distribution system with 
single-phase and three-phase grid-connected photovoltaic (PV) systems using the sequential Monte Carlo technique. A novel 
probabilistic model of output power of PV system is proposed considering random change of weather and hourly variation of 
relevant parameters. The time-varying load with randomness is also incorporated in PLF. Some cases in modified IEEE 33-node 
distribution system have been implemented to verify the method. Time-sequence impacts of PV generation on power flow and 
voltage profiles in distribution system, taking into account the correlation of weather between different zones, are gained from 
results. 
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Solar declination angle. 

Solar altitude angle. 

Solar incident angle. 

Date serial number from Jan. 1st. 

Latitude of PV installation site. 

Longitude of PV installation site. 

Longitude of where standard time corresponds to. 

Hour angle. 

Local time, i.e. true solar time. 

Standard time, refer to Beijing time in this paper. 

Time difference. 

Solar air mass. 

Coefficient of atmospheric transparency. 

Average annual coefficient of atmospheric 

transparency. 

Solar irradiance of exoatmospheric tangential 
plane. 

Solar constant, defined as 1382W/m 2 generally. 

Correction factor of eccentricity from earth to sun. 

Solar direct irradiance through the atmosphere. 

Solar diffuse irradiance through the atmosphere. 

Sunrise time and sunset time. 

Solar irradiance on the inclined PV arrays. 

Inclination of PV arrays to the horizontal plane. 

Length and width of PV arrays. 

Distance between the front and back PV arrays. 

Number of rows of PV arrays. 

Reflectance of the ground. 

Temperature of PV cells. 

Air temperature. 



r\ Efficiency of PV cells. 

r\„ Efficiency of PV cells at reference temperature. 

T r Reference temperature, 25 °C or 298K. 

S Temperature coefficient of PV arrays. 

P s Ratio of shadow. 

S Effective sun light area of PV arrays, 

//inv Efficiency of inverter. 

r/ mpp t Efficiency of maximum power point tracking 
(MPPT). 

Ppv Output power of grid-connected PV system. 

P N Rated power of grid-connected PV system. 

U,U 1,1/2 Uniformly distributed random number. 

T m „ Run time of PV system. 

Repair time of P V system. 

, jU pv Fault rate and repair rate of PV system. 
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1. Introduction 

There was an increased utilization of photovoltaic (PV) 
during recent years due to liberalized markets and the global 
trend of reducing greenhouse gas emissions [1]. The 
development of PV is promoted by favorable policies all over 
the world. Chinese governments are also making efforts to 
encourage the progress of PV. With the development of PV 
generation technique, more and more PV plants and 
family-hold single-phase PV systems are connected to 
distribution networks in China. However, integration of PV 
generation brings various impacts on load flow, voltage 
profile, voltage regulation, short-circuit current, protections, 
et al. in distribution system, since actual distribution networks 
are designed as radially operated and PV generation changes 
the single-source network to multi-source one [2-4]. 
Regarding load flow, the grid-connected PV systems change 



its level, direction and distribution. The uncertainty also 
increases in load flow because the output power of PV 
systems are random due to variation of solar irradiance and 
weather. Based on the change of load flow, sometimes the 
voltage may cross normal limit which means upper limit 
generally. If there are single-phase PV systems connected to 
the network, voltage imbalance which is measured by 
negative sequence voltage unbalance factor (NSVUF) will 
occur and the NSVUF may be beyond stated limit sometimes. 
It is unable to evaluate the impact of stochastic change on the 
operation of distribution system roundly by simply using 
deterministic load flow analysis which aims at specific input 
parameters of the network, but the problem is solved via 
probabilistic load flow (PLF). 

Currently, studies of PLF in distribution system with PV 
generation focus on voltage quality and include two main 
parts, i.e., the output model of grid-connected PV system and 
the analysis method of PLF. Some researches assume that the 
probability density function (PDF) of the output of PV system 
follows Beta distribution [5] or normal distribution [6] as well 
as the solar irradiance. In [7] parameters of PV cells and 
historical weather are used to calculate the output of PV 
system at 15-minute time step and the sequential output 
power is obtained. The random nature of hourly clearness 
index is taken into account and a complex PDF of the output 
of PV system is yielded in [1] [8] [9], while the random 
nature of hourly diffuse index is considered as well as hourly 
clearness index and the PDF of the PV output is gained 
according to the PDFs of the overall solar irradiance and the 
diffuse irradiance in [10]. In addition, dynamic random 
variable is defined in [11] and the probabilistic model of PV 
output is obtained using superposition of the certain rule and 
random variable. As far as the analysis method, the PLF can 
be solved numerically or analytically, so the method of 
combined cumulants and Gram-Charlier [5] [10] [11] or 
Cornish-Fisher [9], two-point estimate method [6] and Monte 
Carlo techniques [1] [7] are used. But PDFs of voltage, mean 
values of voltage, et al. obtained from these studies are static 
without time-sequence change except the work in [11]. 
Moreover, all studies concentrate on three-phase symmetrical 
distribution system with three-phase grid-connected PV 
systems, which are not appropriate for unsymmetrical 
network with single-phase grid-connected PV systems. For 
example, the voltage imbalance caused by single-phase PV 
systems can not be analyzed in above studies. 

Then, to solve these problems, this paper presents 
three-phase PLF for distribution system with three-phase and 
single-phase grid-connected PV systems, using the sequential 
Monte Carlo simulation. In the procedure, the sequential 
probabilistic model of PV output power is illustrated in detail 
considering solar irradiance, random change of weather, 
efficiency and effective sunlight area of PV arrays, efficiency 
of MPPT and inverter, and the reliability of overall 
grid-connected PV systems. Besides, the time-varying model 
of random loads is also taken into account. The method can 
be used to provide power flows and voltage profiles at any 
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time during the 24-hour cycle of a specific day and month, 
and it will be verified to be feasible in modified IEEE 
33-node distribution system with three-phase unsymmetrical 
parameters. 

2. Proposed Probabilistic Model of Pv Output 
Power 

The output power of grid-connected PV system relates to the 
solar irradiance, the conversion efficiency of PV cells, the 
effective sun light area of PV arrays, the efficiency of MPPT 
and inverter, et al. 

2.1. Model of Solar Irradiance 

Solar irradiance includes direct irradiance and diffuse 
irradiance and changes with different places, seasons and 
hours. Therefore, solar irradiance is variable for 24 hours in 
one day, which can be calculated as follows. 

2.1.1 Model of solar irradiance under clear sky 

The model of solar irradiance under clear sky is regular 
without the random effect. The declination, altitude angle, 
solar air mass, coefficient of atmospheric transparency and 
solar irradiance of exoatmospheric tangential plane should be 
calculated. 

According to the Cooper equation, the declination angle is 
given by 

£ = 23.45sin(360- 284 + " ) (1) 
365 

The altitude angle, which is in the form of Sin, is obtained 

by 

sin(«) = sin(^) s'm(S) + cos(^) cos(£> ) cos(<») (2) 
Where, 

<y = 15° -(t-12) (3) 

Where, 

t = t sl +E/60-l5x(L st -L loc ) (4) 

Where, L st is 120° , i.e. longitude of Beijing, in China. 
And E is given by 

E = 9.87 sin(25) - 7.53cos(5) - 1 .5 sin(fl) (5) 

Where, 

„ 3600-81) 

B = W7~, ^ 

364 

According to Kasten experienced equation [12], the solar 
air mass is given by 

1 

m = r^r (7) 

sin(a) + 0.15(a + 3.885r L253 

With the method of MO-1 [13], the coefficient of 
atmospheric transparency is obtained by 



m 



-(lgf, +0.009)/1.547 



(8) 



The solar irradiance of exoatmospheric tangential plane is 
given by 



Where, 



G o = /^ Ic sin(«) 



2,7777 

^ = l + 0.033cos(— ) 
365 



(9) 



(10) 



According to Eq.l-Eq.10, the solar direct irradiance and 
diffuse irradiance through the atmosphere can be obtained by 
Eq.ll and Eq.12 respectively. 



G n =G P m =yG sc sin(a)P" 



1 \-P m 

G A = —?G, r sin(a) 

d 2 ' 5C l-1.41n/> 



(11) 



(12) 



With focus on impacts of grid-connected PV systems on 
voltage profiles, this paper assumes that voltage would not 
cross limit when at night. Thus, sunrise time and sunset time 
should be calculated as follows: 

t rise = 12- arccos(- tan 8 tan q>) 1 1 5 (13) 



t set = 12 + arccos(-tan£tan^)/15 
2.1.2 Solar irradiance on the inclined PV arrays 



(14) 



The solar irradiance on the inclined PV arrays concerns with 
weather, pollutant in the air, dust on the arrays, et al. The 
effect of weather will be elaborated in the next subsection. 
The effect of pollutant is so little that it is ignored. The dust 
could be washed away by rain and self surface cleaning solar 
panels appear presently, so dust is also out of consideration. 
Then, solar irradiance on the inclined PV arrays is given by 



l + co8(^) \J 1-cosQg) 

G = G n R b + G d ( = ) + [ G n + G d )P Z 

L \ L 



Where, 



R b = cos(#)/sin(«) 



(15) 



(16) 



Where, due to the PV arrays usually fixed due south in 
China, we have 

cos(6>) = sin( 5) sin(^> - /?) + cos((p - (3) cos(£) cos(<y) 

(17) 

2.1.3 The effect of random change of weather 

The random change of weather is caused by clouds in reality. 
The clouds may develop, dissipate or vary from one kind to 
another due to the change of meteorology. Different kinds of 
clouds have different impacts on PV systems. Since it is hard 
to get the information of variable clouds, the weather 
becomes the study point instead of clouds. According to the 
analysis of probability statistics, the variation of weather 
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follows some rule that every transformation of state is based 
on the last state and independent of previous states, which 
could be taken as Markov stochastic process [14]. Thus we 
could use Markov transition probability matrix to describe it 
[15]. According to the variation of clouds, the weather could 
be divided into sun day, cloudy day, overcast sky, rain and 
snow for simplified analysis, and the Markov state transition 
diagram of them is presented in Figure 1. 

^22 A, 33 




Figure 1. State transition diagram of weather 

The transition probability (times in one month) between 
two states can be achieved according to historical data of 
weather. Assuming that the present state of weather is i, it 
needs to determine the next state and the transition time, 
which can be done as follows. Giving the present time tj,the 
state transition probability between present weather and the 
five states is Xij, then, the transition time ty is obtained by 



t. 



base 



I; 



U (j= 1,2,3,4,5) 



(18) 



Where, j is the next probabilistic state; t base is the number 
of hours of each month (e.g. months of 30 days are 720 hours; 
February is 672 hours, and left are 744 hours). 

The least value t^,, = Min (ty) is the desired transition time 
and the corresponding j is the next transition state. In this way, 
the sustained time sequence of different weather can be 
obtained and hourly state can be acquired via sequential 
sampling. Next, supposing solar irradiance is unit on sun day, 
the attenuation coefficient (Q of solar irradiance in different 
weather for each hour can be obtained based on the analysis 
of historical data of solar irradiance. At last, the solar 
irradiance on PV arrays for each hour is given by 



(19) 



2.2. The Efficiency of PV Cells 

The conversion efficiency of PV cells is influenced by the 



temperature of PV cells. The change of temperature involves 
wind speed and wind direction, but it is hard to quantify the 
relationships. Thus, the wind is assumed to be irrelevant in 
this paper and the temperature of PV cells is gained by hourly 
air temperature and solar irradiance. 

The hourly air temperature is determined by weather. The 
hourly state of weather in every day is decided by weather 
sampling and the whole day's weather is the state that lasts 
the longest time on daytime. Next, it needs to decide the day's 
highest temperature T max and lowest temperature T^ 
according to historical data of typical weather and then the 
temperature for each hour is obtained by 



Ut) = 



T +T ■ +(T 

max mm V m; 



24 



(20) 



Where, t p is the time when the highest temperature occurs, 
taking it as 14:00 generally. 

The temperature of PV cells is given by 

r = r a + o.03G s (2i) 

Then the efficiency of P V cells is given by 

T] = T] [l-£(T-T r )] (22) 

2.3. The Effective Sun Light Area of PV Arrays 

The back arrays, except the first row, may be affected by the 
shadow of the front arrays with the motion of sun, so the 
effective sun light area should be calculated. 

As shown in Figure 2, x is the length of shadow and is 
obtained by following expressions. 



X = 



H = Lsin(/?) 
d = d - Lcos(fi) 
H - d tan(«) 



tan(«)cos(/?) + sin(/?) 
The ratio of shadow is given by 

s L 

Then the effective sun light area can be obtained by 
S =WxL + WxLx(R-l)x(l-P s ) 

Solar Beam 



(23) 
(24) 

(25) 

(26) 
(27) 
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Figure 2. Schematic diagram of photovoltaic arrays shadow 
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2.4. The Output Power of Grid-connected PV 
System 

Considering the efficiency of MPPT and inverter, the output 
power of PV system is given by 



pv 



(28) 



For grid-connected PV system, when the output power is 
lower than q% of rated power, the system should be separated 
from distribution networks. It keeps to the criterion that when 
P PV is bigger than P N -q%, the PV system operates normally 
and the output is P PV ; otherwise, the PV system closes and the 
output is zero. 

2.5. Reliability Sampling of Grid-connected PV 
System 

In this paper, the reliability of PV system is also considered in 
simulation. Taking the PV system as a whole when studying 
the reliability parameters (i.e. fault rate and repair rate), 
because it is not easy to get the reliability parameters of PV 
arrays, inverter and filter respectively [16]. In simulation the 
state of whole grid-connected PV system is gained by 
sequential sampling. The run time and repair time are given 
by following equations. 



T -*™.V, 

run n ^ 1 

T =*™.U 

-'repair 



(29) 



(30) 



ft 



pv 



When the system is in fault, the output is zero, else the 
output is P PV . 

3. Three-Phase Probabilistic Load Flow in 
Distribution System with Grid- Connected Pv 
Systems 

3.1. Load Model 

The time-varying load model [17] is used and the randomness 
of loads is considered in this paper. The load P L of each hour 
is assumed to follow normal distribution [1]. The PDF of P L 
is given by the following expression: 



HP L ) = 



l 



, 2a 1 



(31) 



Where, p is mean value; <J is standard deviation. 
3.2. Three-phase Load Flow 

The forward/backward substitution method [18] is used to 
solve three-phase load flow. It includes following steps: 

1) Read data of distribution system. 

2) Initialize voltage in the network. 

3) Calculate injected current at nodes. 

4) Sum up line section current via forward sweep. 
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5) Update voltage via backward sweep. 

6) Repeat iteration of 3) to 5) until the convergence criteria 
are satisfied. 

7) Finish the calculation. 

Here, assuming that the grid-connected PV systems operate 
with constant power factor (PF), so they are considered as 
PQ-type nodes in load flow analysis. 

3.3. Probabilistic Load Flow 

The authors developed a program to solve three-phase PLF in 
distribution system with three-phase and single-phase 
grid-connected PV systems based on the simulation method 
presented above, with the help of VC++ software. The loads 
and PV generator productions are assumed to be mutually 
independent. Variations of load flow and voltage profiles 
within a day can be achieved by year or month or day 
simulation in this program. Taking the day type for example, 
the main program flow chart is shown in Figure 3. 



Start 



Read Data 



Input Simulation Date 



Set Times N of Simulation 



Calculate sunrise time t rise 
and sunset time La 



i = i+ 1 
x 



Deal with 
Results 




Deal with 
Probability Results 

( End A ; 



t = t+ 1 



Calculate Output of PV 

i 

PV Systems 
Reliability Sampling 

Calculate Load 

i 

Solve Three-phase Load 
Flow 



Figure 3. Main program flow chart of day type 

4. Case Study 

The IEEE 33-node distribution system, as shown in Figure 4, 
is analyzed in this section, given the terminal voltage 10.5kV. 
Three phase parameters of the network are modified to be 
unsymmetrical, which are introduced in [18]. 




Figure 4. IEEE 33-node distribution system 

4.1. Explanation For Cases 

1) This paper assumes that PV installation site is Hefei 
(Anhui Province, China), whose latitude is 31.86 N and 
longitude is 117.27 E, , taking the day type for example. The 
simulation date is May 30th, on which the sunrise time is 
5:07 and sunset time is 19:03. The time have been 
transformed to Beijing time. 

2) Weather data needed are from local meteorological 
department. The PF of grid-connected PV systems are unit 
and other parameters for PV systems are not given in this 
paper due to space limit. 

3) The power of single-phase PV system at lOkV node is 
assumed to be the sum of the power of each single-phase PV 
system owned by low-voltage consumers that are connected 
to this lOkV node through transformers. 

4) The normal voltage of single phase in lOkV distribution 
system is limited between 9kV and 10.7kV in China. When 
outside the limit, we take it as voltage beyond limit. 

5) Voltage imbalance is discussed in cases. Based on 
Chinese GB/T15543, the NSVUF at grid-connected node 
should be less than 2% for three-phase PV systems. In cases, 
the NSVUF is calculated by 

V 



VUF 



x 100% 



(32) 



Where, V and V + are amplitudes of negative sequence 
voltage and positive sequence voltage respectively. 

6) Weather difference in different zones is considered in 
cases. As shown in Figure 4, if zone A is near zone B, the 
weather in them changes synchronously, thus output power of 
PV systems in them change with the same step, i.e., they are 
correlated. Otherwise, the output powers of PV systems are 
uncorrelated in zone A and B. 

4.2. Cases Design and Result Analysis 

Different cases are designed in order to verify the feasibility 
of the proposed method and study the impact of 
grid-connected PV systems on load flow in distribution 
system. The main points are the impacts of different types of 
grid-connected PV systems as well as the correlation of 
weather on voltage, since voltage quality is primarily 



significant in distribution system. The number of simulations 
is 10000 in order to gain accurate result in each case. 

Case 1 : no grid-connected PV system, which is the base of 
the following comparison. 

When the weather is perfectly correlated in zone A and B, 

Case 2: there are there-phase PV systems rated 300kW at 
node 16, 17, 18, 31, 32, 33; 

Case 3: there is another single-phase PV system rated 
lOOkW connected to phase A at above six nodes respectively 
based on case 2; 

Case 4: there is another single-phase PV system rated 
lOOkW connected to phase B at above six nodes respectively 
based on case 3. 

When the weather is uncorrelated in zone A and B, repeat 
case 2, case 3 and case 4, which are defined as case 5, case 6 
and case 7 respectively. 

The analysis of results just reported power flow of phase A 
in line 17-18 and voltage profiles of node 18 for the sake of 
brevity, since other lines or nodes can be analyzed in a similar 
way. 

The variations of mean value of active power and reactive 
power of phase A in line 17-18 in easel and case2 are 
presented in Figure 5. It can be noted that after integration of 
PV systems the direction of active power is changed from 
7:00 to 17:00 and the amplitude of reverse flow increase first 
up to highest at noon and then decrease, but the direction is 
normal in the early morning and at dusk. Meanwhile, the 
reactive power keeps constant nearly whether PV systems are 
connected to the network or not, because the PF of PV 
systems are unit, i.e., the output reactive power of PV 
systems into the grid is zero. 
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Case2 (A) 
Case2 (R) 
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Figure 5. Variation of mean value of active power and reactive 
power of phase A in line 17-18. (A: active power; R: reactive power) 

The variations of average voltage of each phase at node 1 8 
are shown in Figure 6, Figure 7 and Figure 8, in which 
10.7kV is the upper limit of voltage. They indicate following 
points. Firstly, node voltages are improved with 
grid-connected PV systems and single-phase PV systems 
would improve the voltage further of phase to which the 
single-phase PV systems are connected. Secondly, node 
voltages become higher first and then decrease with 
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integration of PV, and the highest values occur at noon. 
Thirdly, under the same grid-connected condition, when the 
weather is correlated the average voltage is lower than that 
when uncorrelated at the same time, which is obvious around 
noon. 
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Figure 6. Variation of average voltage of phase A at node 1 8 
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Figure 7. Variation of average voltage of phase B at node 18 
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Figure 8. Variation of average voltage of phase C at node 18 



Case2 
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Case4 
Case5 
Case6 
Case? 




9 10 11 12 13 14 15 16 17 18 19 
Hour 

Figure 9. Variation of average NSVUF at node 18 

The variation of average NSVUF at node 18 is presented in 
Figure 9. It indicates points as follows. Firstly, the NSVUF is 
small and varies little when there are only three-phase PV 
systems connected to the network on daytime, but it becomes 
bigger and increases first and then decreases when 
single-phase PV systems are added. Secondly, keeping the 
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capacities of single-phase PV systems constant, the NSVUF 
is bigger when there is only one phase with single-phase PV 
systems than that when there are two phases with 
single-phase PV systems. Thirdly, under the same 
grid-connected condition, when the weather is correlated the 
NSVUF is smaller than that when uncorrected at the same 
time. 

Table.l presents the probabilisties of each phase voltage 
beyond limit at node 18, in which the time when no voltage is 
beyond limit are removed. It shows following points. Firstly, 
the probabilities of each phase voltage beyond limit increase 
first and then decrease, and are up to highest at noon. 
Secondly, the probabilities beyond limit of phase with 
single-phase grid-connected PV systems are bigger than that 
of other phases. Thirdly, under the same grid-connected 
condition, the probabilities of each phase voltage beyond 
limit are smaller when the weather is correlated than that 
when uncorrected at the same time, but the result is opposite 
if no single-phase grid-connected PV system connected to the 
network when the solar irradiation is weak, such as 8:00 and 
16:00. 



Table 1. Variation of probabilities of voltage beyond limit of each phase at node 18 
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0.5364 
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Case3 


A 


0.3694 


0.6412 


0.7042 


0.7542 
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B 
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0.5918 


0.7038 


0.7542 
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Case4 


A 


0.5502 


0.5490 
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0.7976 
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0.5480 


0.5484 
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Case7 
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5. Conclusions 

A method has been proposed for three-phase PLF in 
distribution system with grid-connected PV systems using the 
technique of sequential Monte Carlo simulation, considering 
the stochastic nature of PV generation and loads. In the case 
of PV generation, probabilistic model of output power of 
grid-connected PV system is presented in detail. Based on the 
method, a relevant program is developed. Some cases have 
been implemented using the program to study the impacts of 
different types of grid-connected PV systems and weather 
correlation on the voltage. The conclusions from results of 
cases can be drawn as follows: 

1) The proposed three-phase PLF method is feasible in 
unsymmetrical distribution system with different types of PV 
systems and it can express the time-sequence change of load 
flow due to hourly probabilistic output power of 
grid-connected PV systems. 

2) PV generation may change the direction of power flow 
in distribution system, which would give rise to problems of 
control and operation. 

3) Different types of grid-connected PV systems have 
different impacts on the node voltages of each phase in 
distribution system. It had better to avoid that there is too big 
difference in capacities of single-phase PV systems connected 
to each phase in order to decrease the NSVUF. 

4) It is also suggested that PV systems with high 
penetration should be connected to the distribution system in 
the same zone, which is helpful to decrease the probability of 
voltage beyond limit as well as the NSVUF. 
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